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In this article, locahzed hnear defect modes due to bandgap guidance in two-dimensional pho- 
tonic lattices with localized or non- localized defects are investigated theoretically. First, when the 
defect is localized and weak, eigenvalues of defect modes bifurcated from edges of Bloch bands are 
derived analytically. It is shown that in an attractive (repulsive) defect, defect modes bifurcate out 
from Bloch-band edges with normal (anomalous) diffraction coefficients. Furthermore, distances be- 
tween defect-mode eigenvalues and Bloch-band edges are exponentially small functions of the defect 
strength, which is very different from the one-dimensional case where such distances are quadrat- 
ically small with the defect strength. It is also found that some defect- mode branches bifurcate 
not from Bloch-band edges, but from quasi-edge points within Bloch bands, which is very unusual. 
Second, when the defect is localized but strong, defect modes are determined numerically. It is 
shown that both the repulsive and attractive defects can support various types of defect modes such 
as fundamental, dipole, quadrupole, and vortex modes. These modes reside in various bandgaps 
of the photonic lattice. As the defect strength increases, defect modes move from lower bandgaps 
to higher ones when the defect is repulsive, but remain within each bandgap when the defect is 
attractive, similar to the one-dimensional case. The same phenomena are observed when the defect 
is held fixed while the applied dc field (which controls the lattice potential) increases. Lastly, if 
the defect is non- localized (i.e. it persists at large distances in the lattice), it is shown that defect 
modes can be embedded inside the continuous spectrum, and they can bifurcate out from edges of 
the continuous spectrum algebraically rather than exponentially. 

PACS numbers: 42.70.Qs, 42.65.Tg 



I. INTRODUCTION 



Study of light propagation in periodic media such as photonic crystals and optically induced photonic lattices is of 
great interest to both fundamental physics and applications [1, 2, 3, 4]. It is well known that the unique feature of 
such periodic systems is the existence of bandgap structures in its linear spectrum. Inside the bands, eigenmodes of 
its spectrum are Bloch waves, while inside bandgaps, wave propagation is forbidden because of the repeated Bragg 
reflections. To guide light in periodic media, one of the convenient ways is to introduce a defect into the medium. 
The defect can support linear localized modes (called defect modes) inside bandgaps of the periodic medium. Such 
modes could only propagate along the defect direction (their propagation along other directions is forbidden for being 
inside bandgaps), thus defect guidance (or bandgap guidance) of light is realized. 

Defects and the corresponding defect modes (DMs) have been widely investigated in the field of photonic crystals 
P, 0], where bandgaps and defect modes are in the temporal frequency domains. Recently, reconfigurable optically 
induced photonic lattices in photorefractive crystals with and without defects were successfully generated [H, @, 0, S] • 
In photonic lattices, bandgaps are in the spatial frequency domains. Linear DMs in one-dimensional (ID) photonic 
lattices have been analyzed both theoretically and experimentally in [8, 9, 10]. Nonlinear defect solitons in ID photonic 
lattices have also been theoretically explored pjj, ll2| (see also a study of DMs and defect solitons in ID waveguide 
arrays in pJi]). A more interesting subject is DMs in 2D photonic lattices, where richer light-guiding possibilities can 
arise. Recently, we have observed that a localized defect in 2D lattices can guide not only fundamental modes, but also 
higher-order modes with delicate tail structures [7]. However, our theoretical understanding on such new types of 2D 
DM structures is still far from satisfactory. For instance, it is still unclear what types of DMs a 2D defect can possibly 
support, and how localized these DMs can be in defects of various depths. We should mention that in uniformly 
periodic photonic lattices, a number of nonlinear localized modes such as fundamental, dipole, vortex, dipole-array 
and vortex-array solitons have been reported [5,, J^, J^, [13, 0, [21, J2.] • Solitons in Bessel-ring lattices and 
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2D quasi-periodic lattices have been reported as well [H, [25[ . These nonlinear modes may be regarded as linear 
DMs of the defect induced by solitons themselves. Hence studies of linear DMs in defective lattices and nonlinear 
solitons in uniform lattices can stimulate each other. Indeed, many types of linear DMs we will report in this paper 
do resemble certain types of strongly localized solitons in 2D uniform lattices (see [2^, [2l|, Hj^). However, many 
differences exist between solitons in uniform lattices and linear DMs in defective lattices. For instance, low-amplitude 
solitons in uniform lattices are very broad and occupy many lattice sites, thus the "defects" these solitons generate 
are very different from the single-site defects as created in the experiments of [7] and considered in several sections of 
this paper. Another difference is that, in uniform lattices, centers of solitons can be located at two positions (on-site 
and off-site) in the ID case [26| and four positions at or between lattice sites in the 2D case |22|. However, in defective 
lattices, linear DMs can only be centered at the defect site. Due to these differences, previous results on solitons in 
uniform lattices can not be naively copied over to linear DMs. 

In this article, linear localized DMs in 2D optically induced photonic lattices with localized or non-localized defects 
are systematically studied theoretically. For localized defects, we show numerically that both the attractive and 
repulsive defects can support various types of 2D DMs such as fundamental, dipole, quadrupole, and vortex modes. 
These modes reside in various bandgaps of the photonic lattice. When the defect is weak, DMs are further studied 
analytically by asymptotic methods. We show that in a weak repulsive defect, DMs bifurcate out from Bloch-band 
edges with anomalous diffraction coefficients. The situation is opposite for attractive defects. We also show that 
distances between DM eigenvalues and Bloch-band edges are exponentially small functions of the defect strength, 
which is very different from the ID case where such dependence is quadratic. In addition, we find that some DM 
branches bifurcate not from Bloch-band edges, but from quasi-edge points within Bloch bands, which is very unusual. 
As the defect strength increases, DMs move from lower bandgaps to higher ones when the defect is repulsive, but 
remain within each bandgap when the defect is attractive. The same phenomena are observed when the defect is 
held fixed while the applied dc field (which controls the lattice potential) increases. If the defect is non-localized (i.e. 
it does not disappear at large distances in the lattice), we show that DMs exhibit some new features, such as DMs 
can be embedded inside the continuous spectrum, and they can bifurcate out from edges of the continuous spectrum 
algebraically rather than exponentially. 



II. THE THEORETICAL MODEL 



We begin our study by considering one ordinarily polarized 2D square-lattice beam with a localized defect and 
one extraordinarily polarized probe beam with very low intensity launched simultaneously into a photorefractive 
crystal. The two beams are mutually incoherent, and the defective lattice beam is assumed to be uniform along the 
propagation direction (such defective 2D lattices have been created in our earlier experiments [7|). In this situation, 
the dimensionless governing equation for the probe beam is [Bl, [TO, [l^, - 



iU, + U,, + Uyy - -r-T^TTTK^ = 0- (2-1) 



Here U is the slowly- varying amplitude of the probe beam, {x^y) are transverse distances (in units of D/tt with D 
being the lattice spacing), z is the propagation distance (in units of 2kiD'^ /n'^), Eq is the applied dc field (in units of 

Il = lo cos'^ix) cos^{y)[l + eFoix, y)] (2.2) 



is the intensity function of the photonic lattice (normalized by Id^h-, where Id is the dark irradiance of the crystal and 
I}) the background illumination), /q is the peak intensity of the otherwise uniform photonic lattice, Fd(x, y) describes 
the shape of the defect, e controls the strength of the defect, /cq = 27r/Ao is the wave number (Aq is the wavelength), 
hi = /cone. Tie is the unperturbed refractive index, and rss is the electro-optic coefficient of the crystal. Notice that the 
period of the square lattice has been normalized to be tt. For localized defects which we are considering in this section 
as well as sections Hill and HVl Fd(x,7/) is a localized function (the case of non- localized defects will be considered in 
Sec. (V)). In our numerical computations, for simplicity, we choose Fd(x,?/) to be 0, [13 

Fz)(x,^) =exp[-(x2+y2)Vl28], (2.3) 

which describes a single-site defect. When e > 0, the lattice intensity II at the defect site is higher than that at the 
surrounding regions, and a defect like this is called an attractive defect. Otherwise, the defect is called a repulsive 
one. The defected lattice profiles with e = ±0.9 will be displayed later in the text (see Figs. [5]and[7j). Throughout 
this paper (except in section [Vj), we will choose the lattice peak intensity to be /q = 6. 
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FIG. 1: (a) Diffraction relation of the uniform lattice potential (|2.8p with Eq = 15 and /o = 6 in the reduced Brillouin zone 
along the direction oiF ^ X ^ M ^ F. Shaded: first three Bloch bands, (b) Potential (i.e. effective refractive index change) 
()2.8p in the (x,y) plane, (c) The first Brillouin zone of the 2D lattice in the reciprocal lattice space. 



Localized defect modes in Eq. (|2.ip are sought in the form of 

U (x, y; z) = u{x, y) exp(-z/i2;), (2.4) 

where ji is the propagation constant, and u{x^ y) is a localized function in x and y. After substituting the above 
expression into Eq. (|2.ip , a linear eigenvalue equation for u{x^y) is derived: 

E 

Uxx + Uyy + [/i - —-—r^. t\u = 0. (2.5) 

1 + iL[x,y) 

In the rest of this paper, we will comprehensively analyze DMs in this 2D Schrodinger equation with defective lattice 
potentials by both analytical and numerical techniques. 

Before the analysis of DMs, let us first look at the diffraction relation and bandgap structure of Eq. (|2.5p without 
defects (i.e., e = 0). According to the Bloch theorem, eigenfunctions of Eq. (|2.5p can be sought in the form of 

uix, y) = e^^^^+^^^^G(x, y; k^M). M = , ^2), (2.6) 

where /i = /i(/ci, k2) is the diffraction relation, wavenumbers /ci, k2 are in the first Brillouin zone, i.e. —1 < /ci, A:2 < 1, 
and G(x, y] ki^k2) is a periodic function in x and y with the same period tt (in normalized units) as the uniform lattice 
of (|2.2p . By substituting the above Bloch solution into Eq. (|2.5p (with e = 0), we find that the diffraction relation 
/i(/ci, /C2) of our 2D uniform lattice will be obtained by solving the following eigenvalue problem 

[{d^ + ikif + {dy + ik2f + V{x, y)]G{x, y) = -fiG{x, y) (2.7) 

with the uniform periodic potential 

y{^.y) = — r — — (2.8) 

^ ^ 1 + /OCOS2(X)COS2(7/) ^ ^ 

Figure [1] shows the diffraction relation of Eq. (|2.5p along three characteristic high-symmetry directions {T ^ X ^ 
M ^ r) of the irreducible Brillouin zone. It can be seen that there exist three complete gaps which are named the 
semi- infinite, first and second gaps respectively. They correspond to the white areas in Figure [1] from the bottom to 
the top, separated by the shaded Bloch bands. 

Loosely speaking, a lattice with an attractive defect is like a uniform lattice superimposed with a bright beam (soli- 
ton) under self-focusing nonlinearity, while a lattice with a repulsive defect corresponds to that under self-defocusing 
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nonlinearity. It is well known that spatial solitons can be formed if the diffraction of the medium is normal (corre- 
sponding to ii" {k) > in our current notations), and the nonlinearity of the medium is self- focusing; or the diffraction 
is anomalous {iJi"{k) < 0) and the nonlinearity is self-defocusing. From Figure [TJa), we can see that the diffraction 
relation /C2) is normal at the lower edge of every band, and anomalous at the upper edge of every band. Thus it 
can be anticipated that in a weak attractive defect, DMs will bifurcate out from the lower edge of every Bloch band; 
while in a weak repulsive defect, DMs will bifurcate out from the upper edge of every Bloch band. For ID defects, 
this heuristic argument has been fully confirmed in fid]. For 2D defects, it will be validated in this paper as well. Of 
course, such heuristic arguments can not offer us any insight on quantitative behaviors of DM bifurcations. In the 
ID case, it has been shown that eigenvalues of DMs depend on the defect strength e quadratically when the defect is 
weak [10|. For 2D defects, DM eigenvalues turn out to be exponentially small with the defect strength e (see below), 
which is distinctively different from the ID case. In the next section, we will study in detail how DMs bifurcate from 
edges of Bloch bands both analytically and numerically. Another interesting phenomenon which we will demonstrate 
there is that some DM branches do not bifurcate from Bloch-band edges. Rather, they bifurcate from points inside 
Bloch bands. 



III. DEPENDENCE OF DEFECT MODES ON THE STRENGTH e OF LOCALIZED DEFECTS 

In this section, we consider DMs in the model equations (|2.ip - (j2.2p under localized defects, and investigate how 
they depend on the defect strength e. First, we analytically study how DMs bifurcate from edges of Bloch bands 
when the defect is weak. The method we will use is analogous to one used in |29|] for eigenvalue bifurcations in the 
Schrodinger equation with a weak radially symmetric potential. Second, we numerically determine various types of 
DMs under both weak and strong defects of the form (|2.3p by directly solving the 2D eigenvalue problem (|2.5p . The 
numerical method we use is a power-conserving squared-operator iteration method developed in [30] . 

A. Bifurcations of defect modes from Bloch-band edges in weak localized defects 

Consider the general two-dimensional perturbed Hill's equation 

"ipxx + i^yy + + V{x, y)] = e/(x, y)tlj, (3.1) 

where the unperturbed potential function V{x^ y) is periodic along both the x and y directions with periods Ti and 
T2 respectively, the perturbation to the potential, fix^y)^ is a 2D-localized defect function in the {x^y) plane (i.e. 
f{x,y) ^ as {x^y) 00), and e <C 1. Application of our general results on Eq. ()3.ip to the specific case of Eq. 
(|2.5p will be given in the end of this subsection. Note that the 2D-localization assumption on the function f{x^y) is 
important for the analysis below. If /(x, y) is not 2D-localized, then DM bifurcation behaviors could be quite different 
(see section |V|). 
When e = 0, Eq. ()3.ip admits Bloch solutions of the form 

tP{x, y) = Bn{x, y; ki, k2) = e^^^^+^^2^Gn(x, y; ki, A^s), /i = ^nih, A^s), (3.2) 

where /i = /in(^i7 ^2) is the diffraction relation of the n-th Bloch surface, (/ci, /C2) lies in the first Brillouin zone, i.e. 
— tt/Ti < ki < tt/Ti, — 7r/T2 < k2 < 7r/T2, and Gn(^, ^; ^1, ^2) is a periodic function in both x and y with periods 
Ti and T2 respectively. All these Bloch modes {Bn{x, y; /ci, /C2), (A^i, /C2) ^ the Brillouin zone, n = 1, 2, ... } form a 
complete set [31]. In addition, the orthogonality condition between these Bloch modes is 

/oo POO 
/ B^{x,y; ki, k2)Bn{x,y] ki,k2)dxdy = (27r)^ 5{ki - ki)5(k2 - k2)5{m - n). (3.3) 
-00 J —00 

Here the Bloch functions have been normalized by 

/o^' ^^/o^' dy\Gn{x,y; /ci,fe2)p _ 

T1T2 " ' 

S{.) is the ^-function, and the superscript represents complex conjugation. When e 7^ 0, localized eigenfunctions 
(i.e. DMs) can bifurcate out from edges of Bloch bands into band gaps. Asymptotic analysis of these DMs for e <C 1 
will be presented below. 

Consider bifurcations of defect modes from an edge point ji = fic the n-th Bloch diffraction surface. If two 
or more Bloch diffraction surfaces share the same edge point, each diffraction surface would generate its own defect 
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mode upon bifurcation, thus one could treat each diffraction surface separately. After bifurcation, the propagation 
constant ja of the defect mode will enter the band gap adjacent to the diffraction surface. Without loss of generality, 
we assume that the edge point of this diffraction surface is located at a F symmetry point of the Brillouin zone where 
(/ci, /C2) = (0,0). If the bifurcation point is located at other points such as M or X points, the analysis would remain 
the same. 

When e 7^ 0, defect modes can be expanded into Bloch waves as 

ilj{x,y) = ^ / dki / dk2 an{ki,k2)Bn{x,y;ki,k2), (3.4) 

n^l-J-Tv/Ti J-7r/T2 

where an(/ci, k2) is the Bloch- mode coefficient in the expansion. In the remainder of this subsection, unless otherwise 
indicated, integrations for dki and dk2 are always over the first Brillouin zone, i.e. the lower and upper limits for ki 
and /c2 in the integrals are always {— 7r/Ti,7r/Ti} and {— 7r/T2,7r/T2} respectively, thus such lower and upper limits 
will be omitted below. 

When solution expansion (|3.4p is substituted into the left hand side of Eq. (|3.ip . we get 

X] / / (/>n(/^i, k2)Bn{x, y; ki, fe) dkidk2 = ef{x, y)ip{x, y), (3.5) 

n=l "^"^ 

where 0n(^i, ^2) is defined as 

0n(^l, fe) = (y-n{ki, fe) [/i - /in(^l, ^2)] • (3.6) 

Due to our localization assumption on the defect function fix^y)^ the right hand side of Eq. (|3.5p is a 2D-localized 
function, thus its Bloch-expansion coefficient (j)n{ki^ k2) is uniformly bounded in the Brillouin zone for all values of n 
and e. When solution expansion (|3.4p is further substituted into the right hand side of Eq. ()3.5p and orthogonality 
conditions (|3.3p utilized, we find that (j)n{^ii ^2) satisfies the following integral equation 

<t>n{kuk2) = 777^ T [[ ^^^i%H-W«,„(/Ci,/C2,fcl,fc2)dfclrffc2, (3.7) 

where function Wm,n in the kernel is defined as 

/CO pCO 
/ /(x, y)B*{x, y; ki, k2)Bm{x, y; ki,k2)dxdy. (3.8) 
-co J —oo 

Again, since f{x^y) is a 2D-localized function, Wm,n{ki^k2^ki^k2) is also uniformly bounded for all (/ci,/c2) and 
(^i^fe) points in the Brillouin zone. The fact of functions and Wm,n being uniformly bounded is important in 
the following calculations. Otherwise (when fix^y) is not 2D-localized, for instance), the results below would not be 
valid. 

At the edge point /i = /ic, d /in /dki = d/in/dk2 = 0. For simplicity, we also assume that d'^ /in/dkidk2 = at this 
edge point — an assumption which is always satisfied for Eqs. (|2.ip - (|2.3p we considered in Sec. 2 due to symmetries of 
its defective lattice. If d'^ jin / dkidk2 7^ in some other problems, the analysis below only needs minor modifications. 
Under the above assumption, the local diffraction function near a F-symmetry edge point can be expanded as 

I~in{ki,k2) = /ic + 71^1 + 72^2 + ^(^15 ^1^2, kl), (3.9) 

where 



1 d'^iin 

71 - 



2 dki 



72 



(0,0) ^ ^^2 



(3.10) 

(0,0) 



Since ja = jHc is an edge point, which should certainly be a local maximum or minimum point of the n-th diffraction 
surface, clearly 71 and 72 must be of the same sign, i.e. 7172 > 0. The DM eigenvalue can be written as 

/i = /ic + cr/l^, (3.11) 



where cr = ±1, and < /i(e) <C 1 when e <C 1. The dependence of h on e will be determined next. 
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When Eqs. ()3.9p and (|3.1ip are substituted into Eq. ()3.7p , we see that only a single term in the summation with 
index m = n makes 0{(j)n) contribution. In this term, the denominator ji — /in(^i5^2) is very small near the V 
symmetry point (bifurcation point) (/ci, = (0, 0), which makes it O(0n) rather than O(e^n)- The rest of the terms 
in the summation give 0{e(j)rn) contributions, because the denominators ji — /im(^i7^2) in such terms are not small 
anywhere in the first Brillouin zone. Thus 

0n(^i, k2) = / / ^^%.H-l^n,n(A:i, ^^2, k2)dkdk2 + 0{e^n)- (3.12) 

Since the first term on the right hand side of the above equation is 0{(j)n) rather than 0(e(/)n), it can balance the left 
hand side of that equation. This issue will be made more clear in the calculations below. In order for the denominator 
in the integral of Eq. (|3.12p not to vanish in the Brillouin zone, we must require that 

cr = -sgn(7i) = -sgn(72). (3.13) 

This simply means that the DM eigenvalue /i must lie inside the band gap as expected. 
Now we substitute expressions (j3.9p and (|3.1ip into Eq. (|3.12p . and can easily find that 



K{ki,k2) 
JJ h^^\-fi\kl^\-f2\\kl 



Mki.k2) = j—^ II , ,l\ ^^,, WnAki.k2M.k2)dkidk2^ 0{ec^n)- (3.14) 



The above equation can be further simplified, up to error O(e0n), as 

^n{kiM) = j^M^.^)WnAkiM.^.^) II ---^dkdk^O{e^n)^ (3.15) 

To calculate the integral in the above equation, we introduce variable scalings: ki = \/|7i | ^i, ^2 = \/\l2 \ k- Then 
Eq. (|3.15p becomes 



K{k- 



.k2)= .^""^ (/>n(Q,Q)^n,n(fel,fe2,Q,Q) // ^d^dk ^ 0{e^n). (3.16) 

(27r)^y^7l7J J J h^^kl+kl 



where the integration is over the scaled Brillouin zone with |^i| < tta/IjiI/Ti and |^2| < 7^^/\J2\/T2. This integration 
region can be replaced by a disk of radius — h'^ in the (^i, k) plane, which causes error of 0(e(/)n) to Eq. (|3.16p . 
Over this disk of radius \J\ — h^^ the integral in Eq. (|3.16p can be easily calculated using polar coordinates to be 
— 27rln/i, thus Eq. (|3.16p becomes 

K(kxM) = -77^^^^(/>n(0,0)iyn,n(/Cl,/C2,0,0)+O(e(/)n). (3.17) 

27ry'7i72 

Lastly, we take /ci = /c2 = in the above equation. In order for it to be consistent, we must have 

= - (3.18) 

When this equation is substituted into (|3.1ip . we finally get a formula for the DM eigenvalue /i as 

/i = /ic + (3.19) 

where a is given by Eq. ()3.13p . 

_ 47rcTV7l7i , . 

^-w„.„(o, 0,0,0)' ^-^-'"^ 

and C is some positive constant. Obviously f5 and e must have the same sign, thus e and crVFn,n(0, 0, 0, 0) must 
have the same sign. Since we have shown that a and 71, 72 have opposite signs, we conclude that the condition for 
defect-mode bifurcations from a F-symmetry edge point is that 

sgn[eiyn,n(0, 0,0,0)] = -sgn(7i) = -sgn(72). (3.21) 
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Under this condition, the DM eigenvalue /i bifurcated from the edge point /ic is given by formula (|3.19p . Its distance 
from the edge point, i.e. ji — fici is exponentially small with the defect strength |e|. This contrasts the ID case where 
such dependence is quadratic [10]. The constant C in formula (|3.19p is much more difficult to calculate analytically, 
thus will not be pursued here. 



If the edge point of DM bifurcations is not at a F symmetry point, the DM bifurcation condition (|3.2ip and the 
DM eigenvalue formula (|3.19p still hold, except that 71,72 and Wn,n in these formulas should now be evaluated at 
the underlying edge point in the Brillouin zone. 

Now we apply the above general results to the special system ()2.2p - (|2.5p we were considering in Sec. 2. When 
e <C 1, this system can be written as 



u. 



■XX + Uyy + [/i + V(x, y)\ u = e/(x, y)u + O(e^), (3.22) 



where V{x,y) is given in Eq. ()2.8p . and 



EqIq cos^ (x) cos^ {y)FD (x, y) 
[1 + Iq cos^(x) cos'^ {y)Y 



In this case, W(0, 0, 0, 0) is always negative, thus DM bifurcation condition (|3.2ip reduces to 

sgn(e) = sgn(7i) = sgn(72). (3.24) 

Thus, in an attractive defect, DMs bifurcate out from normal-diffraction band edges (i.e. lower edges) of Figure [TJa) 
where the diffraction coefficients are positive; in a repulsive defect, DMs bifurcate out from anomalous-diffraction 
band edges (i.e. upper edges) of Figure [Tfa) where the diffraction coefficients are negative. This analytical result is 
in agreement with the heuristic argument in the end of Section 2. 

It should be mentioned that at some band edges, two linearly independent Bloch modes exist. Due to the symmetry 
of the lattice in Eq. (|2.2p . these two Bloch modes are related as B{x^ y) and 5(^, x), where B{x^ y) ^ B{y^ x). Under 
weak defects (|2.3p . two different DMs of forms u{x^y) and u{y,x) but with identical eigenvalues fi will bifurcate out 
from these two Bloch modes of edge points respectively. This can be plainly seen from the above analysis. Since these 
DMs have the same propagation constant, their linear superposition remains a defect mode due to the linear nature 
of Eq. (|2.5p . Such superpositions can create more interesting DM patterns. This issue will be explored in more detail 
in the next subsection. 



B. Numerical results of defect modes in localized defects of various depths 

The analytical results of the previous subsection hold under weak localized defects. If the defect becomes strong (i.e. 
|e| not small), the DM formula (|3.19p will be invalid. In this subsection, we determine DMs in Eq. (|2.5p numerically 
for both weak and strong localized defects, and present various types of DM solutions. In addition, in the case of weak 
defects, we will compare numerical results with analytical ones in the previous subsection, and show that they fully 
agree with each other. The numerical method we will use is a power-conserving squared-operator method developed 
in 0. 

In these numerical computations, we fix = 15, Iq =6, and vary the defect strength parameter e from —1 to 1. 
We found a number of DM branches which are plotted in Figure [2] (solid lines). First, we look at these numerical 
results under weak defects (|e| <C 1). In this case, we see that DMs bifurcate out from edges of every Bloch band. 
When the defect is attractive (e > 0), the bifurcation occurs at the left edge of each Bloch band, while when the 
defect is repulsive (e < 0), the bifurcation occurs at the right edge of each Bloch band. Recall that the left and 
right edges of Bloch bands in Figure [2] correspond to the lower and upper edges of Bloch bands in Figure [TJa) (where 
diffractions are normal and anomalous respectively), thus these DM bifurcation results qualitatively agree with the 
theoretical analysis in the previous subsection. We can further make quantitative comparisons between numerical 
values of /i and the theoretical formula (|3.19p . This will be done in two different ways. One way is that we have 
carefully examined the numerical /i values for |e| <C 1, and found that they are indeed well described by functions of 
the form (|3.19p . Data fitting of numerical values of ja into the form (|3.19p gives /3 values which are very close to the 
theoretical values of ()3.2Qp . For instance, for the DM branch in the semi-infinite bandgap in Figure O numerical data 
fitting for < e < 1 gives pnum = 0.1289 {Cnum = 0.4870). The theoretical p value, obtained from Eqs. (|3?2Q|) and 
()3.23p . is Panai = 0.1297, which agrees with Pnum very well. Another way of quantitative comparison we have done 
is to plot the theoretical formula (|3.19p alongside the numerical curves in Figure [2l To do so, we first calculate the 
theoretical value /3 from Eq. ()3.20p at each band edge. Regarding the constant C in the ja formula (|3.19p , we do not 
have a theoretical expression for it. To get around this problem, we fit this single constant from the numerical values 
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FIG. 2: Bifurcations of defect modes with the defect described by Eq. ()2.2p at = 15 and Iq = 6. Sohd hnes: numerical 
results; dashed lines: analytical results. The shaded regions are the Bloch bands. Profiles of DMs at the circled points in this 
figure are displayed in Figs. [3]and|4l 



of /i. The theoretical formulas thus obtained at every Bloch-band edge are plotted in Figure [2] as dashed lines. Good 
quantitative agreement between numerical and analytical values can be seen as well. 

As |e| increases, DM branches move away from band edges (toward the left when e > and toward the right when 
e < 0). Notice that in attractive defects (e > 0), these branches stay inside their respective bandgaps. But in repulsive 
defects (e < 0), DM branches march to edges of higher Bloch bands and then re-appear in higher bandgaps. For 
instance, the DM branch in the first bandgap reaches the edge of the second Bloch band at e —0.70 and re-appears 
in the second bandgap when e ~ —0.76. These behaviors resemble those in the ID case (see Fig. 4 in [10]). One 
difference between the ID and 2D cases seems to be, in ID repulsive defects, when DM branches approach edges of 
higher Bloch bands, the curves become tangential to the vertical band-edge line. This indicates that these ID defect 
modes can not enter Bloch bands as embedded eigenvalues inside the continuous spectrum. This fact is consistent 
with the statement "(discrete) eigenvalues cannot lie in the continuous spectrum" in [321] for locally perturbed periodic 
potentials in ID linear Schrodinger equations. In the present 2D repulsive defects, on the other hand, the curves seem 
to intersect the vertical band-edge lines non-tangentially. This may lead one to suspect that these DM branches might 
enter Bloch bands and become embedded discrete eigenvalues. This suspicion appears to be wrong however for two 
reasons. First, the mathematical work by Kuchment and Vainberg [33] shows that for separable periodic potentials 
perturbed by localized defects, embedded DMs can not exist. In our theoretical model (|2.5p , the periodic potential is 
non-separable, thus their result does not directly apply. But their result strongly suggests that embedded DMs can 
not exist either in our case. Second, our preliminary numerical computations support this non-existence of embedded 
DMs. We find numerically that when a DM branch enters Bloch bands in Figure [21 the DM mode starts to develop 
oscillatory tails in the far field and becomes non- localized (thus not being a DM anymore). 

Now we examine profiles of defect modes on DM branches in Figure [2j For this purpose, we select one representative 
point from each DM branch, mark them by circles, and label them by letters in Figure [2l DM profiles at these marked 
points are displayed in Figs. [3]and[4l The letter labels for these defect modes are identical to those for the marked 
points on DM branches in Figure [2l First, we look at Figure [3l which shows DM profiles in repulsive defects (e < 0). 
We see that DMs at points 'a, b' of Figure [2] are symmetric in both x and with a dominant hump at the defect 
site, and satisfy the relation u{x^y) = u{y^x). These DMs are the simplest in their respective bandgaps, thus we 
will call them fundamental DMs. Notice that these fundamental DMs are sign-indefinite, i.e. they have nodes where 
the intensities are zero, because they do not lie in the semi-infinite bandgap. Although DMs in Figure [3l^a, b) look 
similar, differences between them (mainly in their tail oscillations) do exist due to their residing in different bandgaps. 
These differences are qualitatively similar to the ID case, which has been carefully examined before (see Fig. 8 (b, 
d) in [lOI and Fig. 4 (b, d) in [8]). The DM branch of point 'c' is more interesting. At each point on this branch, 
there are two linearly independent DMs, the reason being that this DM branch bifurcates from the right edge of the 
second Bloch band where two linearly independent Bloch modes exist (this band edge is a X symmetry point, see Fig. 
1). These two DMs are shown in Figure [3l(c, d). One of them is symmetric in x and anti-symmetric in y, while the 
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FIG. 3: Profiles of defect modes in repulsive defects in Figure [2] (a)-(c): DMs at the circled points 'a, b, c' in Figure O with 
(e, = 0.6, 8.33), (—0.8, 10.73), and (—0.6, 12.31), respectively; (d) the coexisting DM of mode (c); (e) and (f): intensity and 
phase of the vortex mode obtained by superimposing modes (c) and (d) with 7r/2 phase delay, i.e. in the form of u(x,y)-\-iu(y, x); 
(g) and (h): intensity and phase of the diagonally-oriented dipole mode obtained by superimposing modes (c) and (d) with no 
phase delay, i.e. in the form of u{x,y) -\-u{y,x). 




FIG. 4: Profiles of defect modes in attractive defects in Figure [2] (i)-(l): DMs at the circled points 'i, j, k, 1' in Figure [21 with 
(e,/i) = (0.8, 4.39), (0.8, 8.12), (0.8, 11.96), and (0.8, 12.57), respectively, (m) and (n): DMs obtained by superimposing mode 
(1) and its coexisting mode with zero and tt phase delays, i.e. in the form of u{x, y) + u{y^ x) and u{x^ y) — u{y, x) respectively. 



other is opposite. They are related as u{x, y) and u{y^ x) with u{x^ y) ^ u{y^ x). These DMs are dipole-like. However, 
these dipoles are largely confined inside the defect site. They are closely related to certain single- B loch- mode solitons 
reported in [22]. 

Since the DM branch of point 'c' admits two linearly independent defect modes, their arbitrary linear superposition 
would remain a defect mode since Eq. (|2.5p is linear. Such linear superpositions could lead to new interesting DM 
structures. For instance, if the two DMs in Figure [3l^c, d) are superimposed with 7r/2 phase delay, i.e. in the form of 
u{x^y) + iu{y,x)^ we get a vortex- type DM whose intensity and phase structures are shown in Figure [Sje, f). This 
vortex DM is strongly confined in the defect site and looks like a familiar vortex-ring. It is qualitatively similar to 



10 



vortex-cell solitons in periodic lattices under defocusing nonlinearity as reported in [22], as well as linear vortex arrays 
in defected lattices as observed in [7]. However, it is quite different from the gap vortex solitons as observed in [20], 
where the vortex beam itself creates an attractive defect with focusing nonlinearity, while the vortex DM here is 
supported by a repulsive defect. If the two DMs in Figure [3l(c, d) are directly superimposed without phase delays, 
i.e. in the form of u{x^y) + u{y^x)^ we get a dipole-type DM whose intensity and phase profiles are shown in Figure 
[3fg, h). This dipole DM is largely confined in the defect site and resembles the DM of Figure [3jc, d), except that its 
orientation is along the diagonal direction instead. This dipole DM is qualitatively similar to dipole-cell solitons in 
uniform lattices under defocusing nonlinearity as reported in [22]. If the two DMs in Figure [Sjc, d) are superimposed 
with TT phase delay, i.e. in the form of u{x,y) — u{y^x)^ the resulting DM would also be dipole-like but aligned 
along the other diagonal axis. Such a DM is structurally the same as the one shown in Figure [3l^g, h) in view of the 
symmetries of our defected lattice. 

Now we examine DMs in attractive defects, which are shown in Figured! At point 'i' in the semi-infinite bandgap, 
the DM is bell-shaped and is strongly confined inside the defect site [see Figure ID^i)]. This mode is guided by the total 
internal reflection mechanism, which is different from the repeated Bragg reflection mechanism for DMs in higher 
bandgaps. On the DM branch of point 'j', there are two linearly independent DMs which are related as u{x,y) and 
u{y,x). One of them is shown in Figure [H^j). This mode is dipole-like and resembles the one in Figure [31(c). Linear 
superpositions of this mode u{x, y) with its coexisting mode u{y^ x) would generate vortex- and dipole-like DMs similar 
to those in Figure [3l^e, f, g, h). In particular, the vortex DM at this point 'j' would be qualitatively similar to the 
gap vortex soliton as observed in [20]. On the DM branch of point 'k', a single DM exists. This is because the left 
edge of the third Bloch band where this branch of DM bifurcates from is located at the M symmetry point of the 
Brillouin zone (see Fig. 1) and admits only a single Bloch mode. The DM at point 'k' is displayed in Figure [D^k). 
This mode is largely confined at the defect site and is quadrupole-like. On the DM branch of point '1', two linearly 
independent DMs exist which are related as u{x^y) and u{y,x). One of them is shown in Figure [31(1). This mode is 
symmetric in both x and y directions and is tripolar-like (with three dominant humps). Unlike DMs at points 'c, j', 
this DM, superimposed with its coexisting mode with 7r/2 phase delay would not generate vortex modes, since this 
DM has non-zero intensity at the center of the defect. However, their superpositions with zero and tt phase delays, 
i.e. in the forms of u{x^ y) + u{y^ x) and u{x^ y) — u{y^ x) would lead to two structurally different defect modes, which 
are shown in Figure [4l[m, n). The DM in Figure [4l[m) has a dominant hump in the center of the defect, surrounded 
by a negative ring, and with weaker satellite humps further out. The DM in Figure [l^n) is quadrupole-like, but 
oriented differently from the quadrupole-like DM in Figure [Ifk). These DMs resemble the spike-array solitons and 
quadrupole-array solitons in uniform lattices under focusing nonlinearity as reported in [22]. 

Most of the DM branches in Figure [2l bifurcate from edges of Bloch bands. Even the branch of point 'b' in the 
second bandgap of Figure [2] can be traced to the DM bifurcation from the right edge of the first Bloch band. But 
there are exceptions. One example is the DM branch of point T in the upper right corner of the second bandgap in 
Figure [2] (we will call it the '1' branch below). This T branch does not bifurcate from any Bloch-band edge. Careful 
examinations show that DMs on this branch closely resemble Bloch modes at the lowest F-symmetry point in the 
third Bloch band (see Figure [1]). Thus this '1' branch should be considered as bifurcating from that lowest F-symmetry 
point when e <C 1. However, we see from Figure [T] that this lowest F-symmetry point is not an edge point of the third 
Bloch band. Even though it is a local edge (minimum) point of diffraction surfaces in the Brillouin zone, it still lies 
inside the third Bloch band. Let us call this type of points "quasi-edge points" of Bloch bands. Then the '1' branch of 
DMs bifurcates from a quasi-edge point, not a true edge point, of a Bloch band! To illustrate this fact, we connected 
the T branch to this quasi-edge F-symmetry point (at e = 0) by dotted lines through a fitted function of the form 
(j3.19p in Figure [2 This dotted line lies inside the third Bloch band. One question we can ask here is: what is the 
nature of this dotted line? Since it is inside the Bloch band, the mathematical results by Kuchment and Vainberg 
[33] suggest that it can not be a branch of embedded DMs. On the other hand, the T DM branch in the second 
bandgap does bifurcate out along this route. Then how does the T branch bifurcate from the quasi-edge point along 
this route? This question awaits further investigation. We note by passing that the F and M points inside the second 
Bloch band are also quasi-edge points (see Figure [Tl[ a)). Additional DM branches may bifurcate from such quasi-edge 
points as well. 

IV. DEPENDENCE OF DEFECT MODES ON THE APPLIED DC FIELD Eo 

In the previous section, we investigated the bifurcation of DMs as the local-defect strength parameter e varies. In 
experiments on photorefractive crystals, the applied dc field Eq can be adjusted in a wide range much more easily 
than the defect strength. So in this section, we investigate how DMs are affected by the value of Eq in Eqs. (|2.ip - (|2.2p 
under localized defects of ()2.3p . When Eq changes, so does the lattice potential. We consider both attractive and 
repulsive defects at fixed values of e = ±0.9 (the reason for these choices of e values will be given below). For other 
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FIG. 5: (a) Profile of the lattice Il^x^v) in Eq. (|22]) with a repulsive localized defect of and e = -0.9 (/o = 6). (b) The 
corresponding DM branches in the (/i, ^o) parameter space. Shaded: Bloch bands. DMs at the circled points are displayed in 
Figure [6l 

fixed e values, the dependence of DMs on Eq is expected to be qualitatively the same. The reader is reminded that 
throughout this section, the /q value is fixed at /q = 6. 

A. The case of repulsive defects 

In 2D photonic lattices with single-site repulsive defects we experimentally created in [7], e —1. At this value of 
6, however, we found that DMs may exist in high bandgaps, but not in low bandgaps (see Figure [2] where Eq = 15). 
Thus for presentational purposes, we take e = —0.9 in our computations below. The corresponding lattice profile of 
Ihix^y) is shown in Figure [51(a), which still closely resembles the defects in our previous experiments [3]. 

For this defect, we found different types of defect modes at various values of Eq. The results are summarized in 
Figs. [Ifb) andlH In Figure [Sjb), the dependence of DM eigenvalues ji on the dc field Eq is displayed. As can be seen, 
when Eq increases from zero, a DM branch appears in the second bandgap at Eq ^ 8.5. When Eq increases further, 
this branch moves toward the third Bloch band and disappears at Eq ^ 15. But when Eq increases to Eq ^ 23.4, 
another DM branch appears in the third bandgap. This branch moves toward the fourth Bloch band and disappears 
at ^0 ^ 36.7. In the fourth gap, there are two DM branches. The upper branch exists when 45.9 < Eq < 64.2, and the 
lower branch exists when 31.6 < < 40.8. As Eq increases even further, DM branches will continue to migrate from 
lower bandgaps to higher ones. These behaviors are qualitatively similar to the ID case (see Fig. 7 in [10]). Some 
examples of DMs marked by circles in Figure [5] are displayed in Figure [6l We can see that DMs at points 'a, b, c' are 
"fundamental" DMs. At point 'd', there are two linearly independent DMs related by u{x^ y) and u{y^ x), one of which 
is displayed in Figure [6](d). These DMs are dipole-like, similar to the 'c' branch in Figure O Linear superpositions 
of these DMs with 7r/2 and phase delays would generate a vortex DM and a diagonally-oriented dipole DM, whose 
amplitude profiles are shown in Figure [6] (e) and (f) (their phase fields are similar to those in Figure [3l(f, h) and thus 
not shown here). 

B. The case of attractive defects 

For attractive defects, we choose e = 0.9 following the above choice of e = —0.9 for repulsive defects. The lattice 
profile at this e value is shown in Figure [71(a). The dependence of DM eigenvalues fi on the dc field Eq is shown 
in Figure [71(b). Unlike DMs in repulsive defects, DMs in this attractive defect exist in all bandgaps (including the 
semi-infinite and first bandgaps). In addition, DM branches here stay in their respective bandgaps as Eq increases, 
contrasting the repulsive case. DM profiles at marked points in Figure [T^b) are displayed in Figure [HI The DM in 
Figure [81(a) is bell-shaped and all positive, and it is the fundamental DM in the semi-infinite bandgap. This DM is 
similar to the one on the 'i' branch of Figure [21 and is guided by the total internal reflection mechanism. The DM in 
the first bandgap is dipole-like and is displayed in Figure [81(b). This DM branch is similar to the 'j' branch of Figure 
[2j On this branch, there is another coexisting DM which is a 90°-rotation of Figure [81(b). Linear superpositions 
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FIG. 6: (a)-(d): Defect modes supported by the repulsive localized defect of Figure [5j a) at the circled points 'a, b, c, d' of 
FigureEIb), with (Eq, fi)^{10.8, 9.61), (27.6, 22.06), (50.0, 37.88), and (35.0, 32.34), respectively, (e) and (f): Intensity profiles 
of the vortex and dipole modes obtained by superposing mode (d) and its coexisting mode with 7r/2 and zero phase delays. 



(b) 




FIG. 7: (a) Profile of the lattice /l(x,2/) in Eq. (|22} with an attractive localized defect of (|23)) and e = 0.9 (/o = 6). (b) The 
corresponding DM branches in the (/x, ^o) parameter space. Shaded: Bloch bands. DMs at the circled points are displayed in 
Figure [8] 

of these two DMs could generate vortex-like and diagonally-oriented dipole-like DMs, as Figure [3l^e, f, g, h) shows. 
In the second bandgap, there are two DM branches: the 'c' branch and the 'd' branch. DMs on the 'c' branch are 
quadrupole-like, and this branch is similar to the 'k' branch of Figure O On this branch, there is a single linearly 
independent defect mode. The 'd' branch is similar to the T branch of Figure [2l On this branch, there are two 
linearly independent DMs which are tripolar-like and orthogonal to each other [see Figure [H^l)]. A superposition of 
these DMs with zero phase delay creates a hump-ring-type structure which is shown in Figure [8l^d) [see also Figure 
m^m)]. A different superposition would produce a quadrupole-type structure as the one shown in Figure [4l(n). In the 
third and fourth bandgaps, there are additional DM branches. These branches exist at Eq values higher than 15, thus 
have no counterparts in Figure [2] (where Eq = 15). Profiles of these DMs are more exotic as can be seen in Figure [8] 
(e) and (f). 
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FIG. 8: Defect modes supported by attractive localized defects of Figure [7j a) at the circled points 'a, b, c, d, e, f of Figure 
Etb), with (^o,/i) = (15, 4.29), (15, 8.00), (15, 11.84), (15, 12.48), (30, 24.23) and (42, 36.69), respectively. 



V. DEFECT MODES IN NON-LOCALIZED DEFECTS 



In previous sections, we exclusively dealt with defect modes in 2D-localized defects. We found that DM eigenvalues 
bifurcate out from Bloch-band edges exponentially with the defect strength e [see Eq. (j3.19p ]. In addition, DMs can 
not be embedded inside Bloch bands (i.e. the continuous spectrum). In this section, we briefly discuss DMs in non- 
localized defects, and point out two significant differences between DMs in localized and non-localized defects. One is 
that in non- localized defects, DM eigenvalues can bifurcate out from edges of the continuous spectrum algebraically, 
not exponentially, with the defect strength e. The other one is that in non-localized defects, DMs can be embedded 
inside the continuous spectrum as embedded eigenmodes. 

For simplicity, we consider the following linear Schrodinger equation with separable non-localized defects. 



Uxx + Uyy + {Vd{x) + Voiy)} U = -IIU, 

where function Vd is a one-dimensional defected periodic potential of the form 



Vd{x) = 



En 



Iocos^x){l + eFD{x)}' 



(5.1) 



(5.2) 



Fd{x) = exp(— x^/128) is a defect function as has been used before and e is the defect strength. We chose this 
separable defect in Eq. ()5.ip because all its eigenmodes (both discrete and continuous) can be obtained analytically 
(see below). The non- localized nature of this separable defect can be visually seen in Figure [9l^a), where the defected 
potential Vd{x) + Vniy) for Iq = S^Eq = 6, and e = —0.7 is displayed. Note that this potential plotted here differs 
from the usual concept of quantum- mechanics potentials by a sign, and it resembles the refractive-index function in 
optics. We see that along the x and y axes in Figure [91(a), the defect extends to infinity (thus non-localized). The 
above Iq and Eq parameters were chosen because they have been used before in the ID defect-mode analysis of [loj, 
and such ID results will be used to construct the eigenmodes of the 2D defect equation (|5.ip below. 

Since the potential in Eq. (|5.ip is separable, the eigenmodes of this equation can be split into the following form: 



u{x, y) = Ua{x)ub{y), ja = jiia ^ l^h, 
where Ua^u^j^ fia, fJ^b satisfy the following one-dimensional eigenvalue equation: 

Ua,xx + VD{x)Ua = -/J^a^a, 



(5.3) 



(5.4) 



Ub,yy + VD{y)Ub = -jJ^bUb- 



(5.5) 
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These two ID equations for Ua and 1x5 are identical to the ID defect equation we studied comprehensively in [lO|. 
Using the splitting ()5.3p and the ID results in [10], we can construct the entire discrete and continuous spectra of Eq. 

m- 

Before we construct the spectra of Eq. (|5.ip . we need to clarify the definitions of discrete and continuous eigenvalues 
of this 2D equation. Here an eigenvalue is called discrete if its eigenfunction is square-integrable (thus localized along 
all directions in the (x^y) plane). Otherwise it is called continuous. Note that an eigenfunction which is localized 
along one direction (say x axis) but non-localized along another direction (say y axis) corresponds to a continuous, 
not discrete, eigenvalue. 

Now we construct the spectra of Eq. ()5.ip for a specific example with £^0 = 6, /q = 3 and e = 0.8. At these 
parameter values, the discrete eigenvalues and continuous-spectrum intervals (ID Bloch bands) of the ID eigenvalue 
problem (|5.4p are (see p^of ) 

{Ai, A2, A3, . . . } = {2.0847, 4.5002, 7.5951, ...,}, (5.6) 

{[/i,/2], [h.hi [/5,/6],...} = {[2.5781, 2.9493], [4.7553, 6.6010], [7.6250, 11.8775],...}. (5.7) 

Using the relation (|5.3p . we find that the discrete eigenvalues and continuous-spectrum intervals of the 2D eigenvalue 
problem (|5.ip are 

{/ii, /i2, /is, /i4, . . . } = {2Ai, Ai + A2, 2A2, Ai + A3, . . . } = {4.1694, 6.5849, 9.0004, 9.6798, ...,}, (5.8) 

and 

{/^continuum} = {[^1 + 2/2], [Ai + 00)} = {[4.6628, 5.8986], [6.8400, 00)}. (5.9) 

Note that at the left edges of the two continuous-spectrum bands (/i = 4.6628 and 6.8400), the eigenfunctions are 
non-localized along one direction but localized along its orthogonal direction, thus they are not the usual 2D Bloch 
modes (which would have been non- localized along all directions). This is why for Eq. ()5.ip . we do not call these 
continuous-spectrum bands as Bloch bands. 

Repeating the same calculations to other e values, we have constructed the whole spectra of Eq. (|5.ip in the (/i, e) 
plane for Iq = 3 and £^0 = 6. The results are displayed in Figure [9](b). Here solid curves show branches of DMs 
(discrete eigenvalues), and shaded regions are the continuous spectrum. This figure shows two significant features 
which are distinctively different from those in localized defects. One is that several DM branches (such as the 'c' and 
'd' branches) are either partially or completely embedded inside the continuous spectrum. This means that in non- 
localized defects, embedded DMs inside the continuous spectrum do exist. This contrasts the localized-defect case, 
where the mathematical results of Kuchment and Vainberg [s^ and our numerics show non-existence of embedded 
DMs. Another feature of Figure [9jb) is on the quantitative behavior of DM bifurcations from edges of the continuous 
spectrum at small values of the defect strength e. We have shown before that in the ID case, DMs bifurcate out from 
Bloch-band edges quadratically with e [10]. For the 2D problem (|5.ip . using the relation ()5.3p . we can readily see that 
when DMs bifurcate out from edges of the continuous spectrum (see the 'a, b, d' branches for instance), the distance 
between DM eigenvalues and the continuum edges also depends on e quadratically. This contrasts the localized-defect 
case, where we have shown in section Hill that DMs bifurcate out from Bloch-band edges exponentially. 

Even though the defect in Eq. (|5.ip is non-localized rather than localized, the DMs it admits actually are quite 
similar to those in localized defects. To demonstrate, we picked four representative points on the DM branches of 
Figure [9^b). These points are marked by circles and labeled by letters 'a, b, c, d' respectively. Profiles of DMs at 
these four points are displayed in Figure [TOl Comparing these DMs with those of localized defects in Figs. [3] and [H 
we easily see that they are quite similar. In particular, the 'a, b, c, d' branches in Figure [9](b) resemble the 'i, j, k, a' 
branches in Figure [21 as DMs on the corresponding branches are much alike. Notice that DMs in the non-localized 
defect case are a little more spread out than their counterparts in the local-defect case (compare Figs. HJ^k) and fTOT c). 
for instance). This is simply because in the non- localized defect, we chose ^0 = 6, while in the localized defect, 
we chose a much higher value of £^0 = 15. Higher £^0 values induce deeper potential variations, which facilitates 
better confinement of DMs. We need to point out that, even though the DMs at points 'c, d' are embedded inside 
the continuous spectrum, their eigenfunctions are perfectly 2D-localized and square-integrable (see Figure [TOlf d) in 
particular). Thus, the embedded nature of a defect mode does not necessarily cause it to spread out more. 

VI. SUMMARY 

In summary, we have thoroughly investigated defect modes in two-dimensional photonic lattices with localized or 
non- localized defects. When the defect is localized and weak, we analytically determined defect-mode eigenvalues 
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FIG. 9: (a) Profile of a 2D potential with non- localized defects Vd{x) + Vd(2/) in Eqs. ()5.1|) - ()5.2p with Eq = 6, /q = 3, and 
e = —0.7. (b) DM branches supported by this non-localized defect in the (/x, e) parameter space. Shaded: the continuous 
spectrum. DMs at the circled points are displayed in Figure [TOl 
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FIG. 10: Defect modes in non- localized defects (|5.ip - ()5.2|) at the circled points 'a, b, c, d' in Figure El^b), with (e,/i) = (0.8, 
4.17), (0.8,6.59), (0.8, 9.00), and (-0.7,8.61) respectively. 

bifurcated from edges of Bloch bands. We found that in an attractive (repulsive) defect, defect modes bifurcate out 
from Bloch-band edges with normal (anomalous) diffraction coefficients. Furthermore, distances between defect-mode 
eigenvalues and Bloch-band edges are exponentially small functions of the defect strength, which is qualitatively 
different from the ID case. Another interesting phenomenon we found was that, some defect-mode branches bifurcate 
not from Bloch-band edges, but from quasi-edge points within Bloch bands. When the defect is localized but strong, 
defect modes were studied numerically. It was found that both the repulsive and attractive defects can guide various 
types of defect modes such as fundamental, dipole, tripolar, quadrupole, and vortex modes. These modes reside in 
various bandgaps of the photonic lattice. As the defect strength increases, defect modes move from lower bandgaps 
to higher ones when the defect is repulsive, but remain within each bandgap when the defect is attractive. The 
same phenomena are observed when the defect is held fixed while the applied dc field increases. If the defect is non- 
localized (i.e. the defect does not disappear at infinity in the lattice), we showed that DMs exhibit two new features: 
(i) DMs can be embedded inside the continuous spectrum; (ii) DMs can bifurcate out from edges of the continuous 
spectrum algebraically rather than exponentially. These theoretical results pave the way for further experimental 
demonstrations of various type of DMs in a 2D lattice. 
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